Task 1: making an interactive map, with tmaps

First I’ll read in the CA counties data

ca_counties <- read_sf(here("ca_counties","CA_Counties_TIGER2016.shp")) %>% 
  rename(county_name = NAME, land_area = ALAND) 
  
ca_counties %>% st_crs()
## Coordinate Reference System:
##   User input: WGS 84 / Pseudo-Mercator 
##   wkt:
## PROJCRS["WGS 84 / Pseudo-Mercator",
##     BASEGEOGCRS["WGS 84",
##         DATUM["World Geodetic System 1984",
##             ELLIPSOID["WGS 84",6378137,298.257223563,
##                 LENGTHUNIT["metre",1]]],
##         PRIMEM["Greenwich",0,
##             ANGLEUNIT["degree",0.0174532925199433]],
##         ID["EPSG",4326]],
##     CONVERSION["Popular Visualisation Pseudo-Mercator",
##         METHOD["Popular Visualisation Pseudo Mercator",
##             ID["EPSG",1024]],
##         PARAMETER["Latitude of natural origin",0,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8801]],
##         PARAMETER["Longitude of natural origin",0,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8802]],
##         PARAMETER["False easting",0,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8806]],
##         PARAMETER["False northing",0,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8807]]],
##     CS[Cartesian,2],
##         AXIS["easting (X)",east,
##             ORDER[1],
##             LENGTHUNIT["metre",1]],
##         AXIS["northing (Y)",north,
##             ORDER[2],
##             LENGTHUNIT["metre",1]],
##     USAGE[
##         SCOPE["unknown"],
##         AREA["World - 85°S to 85°N"],
##         BBOX[-85.06,-180,85.06,180]],
##     ID["EPSG",3857]]
## This is using EPSG 3857

## Reading in Oil Spill shapefile (I think this is what we'll use)
oil_spill <- read_sf(here("oil_spill_data", "Oil_Spill_Incident_Tracking_%5Bds394%5D.shp")) 

oil_spill %>% st_crs()
## Coordinate Reference System:
##   User input: WGS 84 / Pseudo-Mercator 
##   wkt:
## PROJCRS["WGS 84 / Pseudo-Mercator",
##     BASEGEOGCRS["WGS 84",
##         DATUM["World Geodetic System 1984",
##             ELLIPSOID["WGS 84",6378137,298.257223563,
##                 LENGTHUNIT["metre",1]]],
##         PRIMEM["Greenwich",0,
##             ANGLEUNIT["degree",0.0174532925199433]],
##         ID["EPSG",4326]],
##     CONVERSION["Popular Visualisation Pseudo-Mercator",
##         METHOD["Popular Visualisation Pseudo Mercator",
##             ID["EPSG",1024]],
##         PARAMETER["Latitude of natural origin",0,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8801]],
##         PARAMETER["Longitude of natural origin",0,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8802]],
##         PARAMETER["False easting",0,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8806]],
##         PARAMETER["False northing",0,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8807]]],
##     CS[Cartesian,2],
##         AXIS["easting (X)",east,
##             ORDER[1],
##             LENGTHUNIT["metre",1]],
##         AXIS["northing (Y)",north,
##             ORDER[2],
##             LENGTHUNIT["metre",1]],
##     USAGE[
##         SCOPE["unknown"],
##         AREA["World - 85°S to 85°N"],
##         BBOX[-85.06,-180,85.06,180]],
##     ID["EPSG",3857]]
## This date mathces the crs of ca_counties, with EPSG 3857

Exploratory plotting

ggplot(data = ca_counties) +
  geom_sf(aes(fill = land_area), color = "white", size = 0.3) +
  theme_void() +
  scale_fill_gradientn(colors = c("lightgreen", "green","darkgreen"))

I’ll do the same for the spill data

ggplot(data = oil_spill) +
  geom_sf() +
  theme_void()

tmap_mode(mode = "view")
## tmap mode set to interactive viewing
tm_shape(ca_counties) +
  tm_fill("land_area", palette = "Reds") +
  tm_shape(oil_spill) +
  tm_dots(col = "blue3",
          border.alpha = 0.5,
          alpha = 0.3)

Part 2

## Counting oil spills by county
# 
oil_spill <- oil_spill%>% 
  rename(county_name = LOCALECOUN)

oil_spill <- oil_spill %>% 
  rename(county = county_name)

oil_spill_inland <- oil_spill %>%
  filter(INLANDMARI == "Inland")

# oil_count <- oil_spill %>% 
#   count(county_name)

# county_count <- ca_counties %>% 
#   count(county_name)

## Joining the two df, 

spill_by_county <- ca_counties %>% 
  st_join(oil_spill_inland)

joint_spill <- spill_by_county %>% 
  count(county_name)

Now we have a df with just county names, number of oil spills per county, and the sticky geometry for each inland oil spill.

Now I’ll graph is in ggplot, making a chloropleth!

ggplot(data = joint_spill) +
  geom_sf(aes(fill = n), color = "white", size = 0.2) +
  scale_fill_gradientn(colors = c("yellow","tomato","purple4")) +
  theme_void() +
  labs(fill = "Number of Inland Oil Spills",
       title = "2008 Oil Spills Throughout California")